#!/usr/bin/env python 
# -*- coding: utf-8 -*-

from Bio import SeqIO
from optparse import OptionParser
import sys
from Bio.SeqUtils import GC
from Bio.Alphabet import IUPAC
from Bio.Seq import Seq
from Bio.SubsMat import MatrixInfo as matlist
from Bio import pairwise2
from Bio.pairwise2 import format_alignment
from Bio.Align.Applications import ClustalwCommandline
from Bio.SeqRecord import SeqRecord
import subprocess
import sys
from Bio import AlignIO

def types_id(files):
    types = {}
    for elem in files:
        if elem.endswith("fasta"):
            types[elem] = "fasta"
        elif elem.endswith("gbk") or elem.endswith("gb"):
            types[elem] = "genbank"
    return types

def count_GC(files, types):
    dictGC = {}
    for elem in files:
        gc_values = [GC(rec.seq) for rec in SeqIO.parse(elem, types[elem])]
        dictGC[elem] = gc_values
    return dictGC

def alignment(proteins, matrix, outfileName):
    for i in range(len(proteins)):
        for j in range(i):
            listOfAl = list(pairwise2.align.globaldx(proteins[i], proteins[j], matrix))
            with open(outfileName, "a") as f:
                f.write(format_alignment(*listOfAl[0]))
            print(format_alignment(*listOfAl[0]))

def mult_alignment(inputfile, matrix, outfileName):
    cline = ClustalwCommandline("clustalw2", infile="proteins.fasta", seqnos="ON", matrix="BLOSUM", gapopen=2, gapext=0.5, outfile=outfileName + ".aln")
    sys.stdout = cline()
    input_handle = open(outfileName + ".aln", "rU")
    output_handle = open(outfileName + ".phy", "w")
    alignments = AlignIO.parse(input_handle, "clustal")
    AlignIO.write(alignments, output_handle, "phylip")
    output_handle.close()
    input_handle.close()
    

def convert_to_prots(all_RNAs):
    proteins = []
    for elem in all_RNAs:
        prot = elem.translate().split("*")
        for i in range(1,len(prot)):
            prot[0] += prot[i]
        proteins.append (prot[0])

    length = 0
    with open ("proteins.fasta","wb") as f:
        for i in range(len(proteins)):
            tmp_line = ">" + str(i) + "\n"
            f.write(tmp_line)
            length += len(tmp_line)
            tmp_line = str(proteins[i]) + "\n"
            f.write(tmp_line)
            length += len(tmp_line)
    return proteins, length

def convert_to_RNAs(files, types):
    all_RNAs = []
    for elem in files:
        RNAs = [rec.seq.transcribe() for rec in SeqIO.parse(elem, types[elem])]
        all_RNAs += sorted(RNAs)
    return all_RNAs

def main(files):
    types = types_id(files)
    dictGC = count_GC(files, types)
    files_string = ""
    for f in files:
        files_string += f + " "
    print "Analyse of GC content in files: " + files_string
    print dictGC
    print
    proteins = []
    length = 0
    all_RNAs = convert_to_RNAs(files, types)
    print "RNA sequences in files: " + files_string
    for elem in all_RNAs:
        print elem
    print 
    proteins, length = convert_to_prots(all_RNAs)
    print "Proteins that was transcribed from files (without stop-codons): " + files_string
    for elem in proteins:
        print elem
    print
    matrix = matlist.blosum62
    f = open("align_output.aln", 'w')
    f.close()
    # alignment(proteins, matrix, "align_output.aln")
    mult_alignment(proteins, matrix, "align_output")
            

    
if __name__ == "__main__":
    parser = OptionParser(version = "%prog 1.0", conflict_handler="error", description =
                            """This script makes 4 steps with unittesting using BioPython:
                            1) count GC\n
                            2) transcription from DNA to RNA \n
                            3) translation from RNA to the sequence of aminoacids \n
                            4) alignment of sequences using BLOSUM62 \n
                            USAGE: python2 task2.py -f file_1_name -f file_2_name \n""")
    cmdArgs = sys.argv[1:]
    if len(cmdArgs) == 0:
        print parser.description
    else:
        parser.add_option("-f", "--files", action="append")
        opts, args = parser.parse_args()
        main(opts.files)

